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Abstract. A one-dimensional turbulent convection model 
in the form of a time-dependent diffusion equation for 
the turbulent energy is incorporated into our numerical 
pulsation code. The effect of turbulent convection on the 
structural rearrangement of the static equilibrium star is 
taken into account, and the linear eigenvalues (periods and 
growth-rates) for the complete turbulent convective nona- 
diabatic pulsations of the equilibrium models are calcu- 
lated. The linearized code allows us to perform efficiently 
a systematic survey of the stability of Cepheid models for 
a wide range of astrophysical parameters (L,M ,T e ff) and 
the turbulent convective parameters. 

A sensitivity analysis of the properties of several se- 
quences of Galactic Cepheids is performed with the goal 
of using observational constraints, such as the shape and 
location of the fundamental and first overtone instabil- 
ity strips, to calibrate the turbulent convective input pa- 
rameters. The locations of the red and blue edges of the 
fundamental and first overtone instability strips and the 
overall strip widths are largely determined by the param- 
eters that determine the mixing length, and strengths of 
the eddy viscosity and of the convective flux. The remain- 
ing parameters can be used for fine-tuning agreement with 
observations. 
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1. Introduction 

Cepheid studies have acknowledged, right from the begin- 
ning, that convection should occur in Cepheid envelopes, 
and that they should alter the properties of both the pul- 
sations and the underlying equilibrium models. Since the 
Prandtl number characteristic of stellar material is very 



small, Pr ~ O(10~ 9 ), and the Rayleigh numbers associ- 
ated with the partial ionization regions (PIR) are enor- 
mous, Ra <; O(10 12 ), stellar convection is far into the 



hard turbulent regime (Krishnamurti & Howard 1983) 
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The cooler Cepheid stars are believed to have the more 
extensive convection zones, leading to pulsational stability 
beyond the red edge of the standard instability strip (IS). 
Because of the estimated inefficiency of convection and be- 
cause of the difficulties of treating turbulent convection, 
a lot of work has disregarded it altogether, the hope hav- 
ing been that, except for providing a red edge, convection 
would play a subdued role. And purely radiative models 
have indeed provided a reasonable overall agreement with 
observations. 

More recently it has become increasingly clear that an 
impasse has been reached with radiative models, and that 
it is not possible to improve the agreement with observa- 
tions (e.g. Buchler 1997). In particular, the amplitudes of 
pulsation are systematically too large. No consistent set 
of pseudo-viscous parameters can correctly limit the am- 
plitudes of both fundamental and first overtone pulsations 
without adversely affecting the stability of the correspond- 
ing limit cycles, and the agreement with observation. The 
new EROS and MACHO data of the lo w met allicity Mag- 
ellanic Cloud Cepheids (Buchler et al. 1996 ) reveals fur- 
ther disagreement between stellar evolution and radiative 
stellar pulsation models. Clearly some additional physi- 
cal dissipation is missing in the radiative codes (see also 
Kovacs 1990), and turbulence and convection are the pri- 
mary suspects. 

Simulation of even the simplest turbulent convection 
(TC) problem is a formidable numerical challenge; to ad- 
equately resolve the many time- and length-scales of fully 
turbulent 3D convection at stellar Rayleigh and Prandtl 
numbers is totally out of question, and astrophysicists 
have been looking for ID recipes to give an approximate, 
but acceptable treatment of this phenomenon. Thus the 
1925 mixing length picture of Prandtl was readily adapted 
to stellar envelopes by Bohm-Vitense, and subsequently 
reformulated in many variations, all of which are equiva- 
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lent for time-independent problems (for an overview, see 
Baker |1987| ). But time-independent MLT was never able 
to resolve the red edge problem, mostly because it is time- 
dependent dissipation introduced by eddy viscosity that 
provides a clear red edge, as we will show in § 4.6. Spiegel 



1963 first extended MLT to the time-dependent and nonlo- 



cal cases. The earliest attempts to include a time-dependent 



mixing length model in pulsation (Gough 1977, Unno 1967) 



were not successful because they were too local in space. 
Subsequently, Gough's theory was extended to be less lo- 
cal, and eventually was used in a full linearization of Solar 
modes (Balmforth 1992 ), analogous to what we h ave d one 
here. Unno's theory was also developed (Castor 1968) as 
a diffusion model, and simplified to a si ngle eq uation for 
the turbulent energy e t by Stellingwerf ( 1982 ). More re- 
cent applications of this and similar turbule nt di ffusion 
recipes have b een m ade by Gonczi & Osaki 198C, Bono 
& Stellingwerf 1994 , Gehm eyr & Winkler 1992j , hereafter 
GW) and Kuhfuss ( 1986 ). Further closely related mod- 
els also appear in the literature, such as those of Canuto 
( 1991 ) but have yet to be applied. 

All ID model TC equations contain several dimension- 
less, order unity parameters that are directly related to 
the physical quantities of the model. In principle, the va- 
lidity of the ID model equation can be checked against 
experiments or against detailed 3D simulations of turbu- 
lent convection, and values of the unknown parameters 
can be extracted. However, neither such experiments nor 
such simulations are presently available. We are obliged 
to calibrate these parameters indirectly with the help of 
observational astronomical constraints, as suggested more 
than a decade ago by Stellingwerf (1982). In the best out- 
come, observational constraints could select one model in 
favor of the others. On the other hand we might discover 
that single equation models for turbulent convecti on are 
insufficient, or that plumes (e.g. Rieutord & Zahn 1995 ) 
play an essential role in the convective transport. 

This programme has recently been initiated by incor- 
porating these recipes into hydrodynamic pulsation codes. 
Fully nonlinear models have been used to determine the 
stability properties of convective pulsators (e.g. Bono & 
Stellingwerf |1994| , Bono & Marconi [1997D , but this ap- 
proach is inefficient, and unsuitable when a large number 
of models is required. The approach we follow in this pa- 
per is to compute the linear nonadiabatic properties di- 
rectly by linearizing the hydrodynamic equations about 
the equilibrium model and by then computing the linear 
eigenvalues, a procedure which is relatively fast. This al- 
lows us to perform an extensive survey of the sensitivity 
to the TC parameters. 

2. The Turbulent Convective (TC) Model 
Equations 

We have adopted a TC model equation that consists of 
a single diffusion equation for the turbulent energy and 



that includes the effects of the turbulent stress, turbulent 
energy flux and turbulent energy production on the mean 
momentum and energy; our model is essentially that used 
by Gehmeyr (1992). In the following, we reproduce our 
PDE's for the purpose of notation and for completeness. 
In addition to dr/dt = u we have the momentum equation, 
modified to include the effects of turbulent stress, that are 
cast here in the form of a turbulent pressure and an eddy 
pressure 



du A o d , , GM(r) 



(1) 



The pressure p and specific energy e represent the com- 
bined gas and radiation quantities. We use the specific 
internal energy equation in the form 



dT__f (de\ \ d(r 2 u) 
v dt I V dv ) t ) dm 



_d_ 

dm 



[r 2 (F c + F r )] +C. (2) 



The turbulent specific energy itself satisfies a diffusion 
equation of the form: 



dt dm 1 J 

where, 



(Pt 



^ d(r 2 u) 
dm 



C = a d 



1/2 

a\H p 



(e t -S t - e Q ) . 



C, (3) 



(4) 



Following GW we have included a negligible, nonzero back- 
ground turbulent energy, e , which is set equal to a small 
constant (e = 10 4 erg cm -3 ) value to avoid numerical 
problems. The turbulent and eddy pressures are given by: 



Pt = a p p e t 



TT l/2 dU 

Vv = -ci v aAli p pe t ' — 



(5) 



(6) 



The mixing length A appears as the product A = a\H p to 
emphasize that it is taken to be some fraction a\ of the 
pressure scale height H p = pr 2 /(pGM). The convective 
and turbulent fluxes are: 



F c = a c a\ p e\' 2 c p T Y , 



r u 1/2 de t 

r t — -a t a\H p e t . 

The source of turbulence that appears in C is given by 

St = a s a A (e t pf3TY/ P y /2 , (9) 

where the instability criterion has been embedded into the 
convenient, dimcnsionlcss entropy gradient: 



(7) 

(8) 



Y 



H p ds 
c„ dr 



(10) 
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Here (3 is the thermal expansion coefficient, and all other 
symbols have their usual meanings. 

In the absence of turbulent diffusion, the formulations 
of Stellingwerf (1982), of GW, and ours are equivalent; all 
three reduce to time-dependent local MLT in this approx- 
imation. With the diffusion term included, however, the 
formulations are no longer identical. As Eq. shows we 
have taken the convective flux to be directly proportional 
to the entropy gradient ds/dr, in the manner of Gehmeyr; 
but in Eq. (^) we have chosen to describe the turbulent 
energy as in Stellingwerf. While this would represent an 
inconsistency in entirely local MLT, here it simply implies 
a different evaluation of the correlation (w'T') in the def- 
inition of the source term, St- In single equation models 
such as this, the absence of a dynamic equation for (w'T') 
necessitates such ad hoc evaluations. 

The TC recipe contains seven dimensionless parame- 
ters of order unity denoted by a's. It is important to test 
the sensitivity of the results to variations of the a's in some 
reasonable range, say 0.1 & a ^5 10. However, these param- 
eters also serve as direct connections to the physical inputs 
of the model. For example: otd, in the absence of diffusion, 
allows us to adjust the timescale of the turbulent energy 
growth and decay (GW); a a controls the mixing length as 
a function of the local pressure scale height; a„ represents 
a constant eddy viscosity; and a p — (74 — 1) is related to 
the adiabatic compressibility index j t of turbulent eddies. 

In our numerical study Eqs. (0)-(|^) are differenced in 
the spirit of the Fraley scheme (see Stellingwerf 1974 and 
Buchler, Kollath & Marom 1998[ ) using the Lagrangean 
mass as the radial variable. 

3. Static Models 

In Cepheids the PIRs are the source both of convection 
and of pulsational driving. Because convection causes a 
structural rearrangement of the stellar model, convective 
stars can have different properties from their radiative 
counterparts even if the direct interplay between convec- 
tion and pulsation is disregarded. In our linearized code it 
is possible to adjust parameters of either the equilibrium 
model, the linearized model, or both; we can therefore 
show the separate effects of the equilibrium structure and 
of the interaction of convection and pulsation on the sta- 
bility of the star 



4.6 



We illustrate this dichotomy in 
by decoupling the turbulent and hydrodynamic pertur- 
bations. But first, we discuss the properties of the static 
equilibrium models. 

The construction of the TC equilibrium models pro- 
ceeds in two steps: first, we integrate the time-independent 
form of Eqs. ([j]-|]) from the surface inward, but in the ab- 
sence of turbulent diffusion (a t =0). This integration be- 
comes the initial guess from which we solve the complete 
time-independent equation by Newton-Raphson iteration. 
(As an alternative we could iteratively solve Eqs. |l] - ||) 
followed by Eq. |§) 




L |VpJ/|v Plolal | 








-a 

Log(p/p cenle 

Fig. 1. Standard Cepheid model sequence (M=5Mq, 
L=2090L o , a c = 3, a a = 0.375) with T ef f(m K) = 
4300, 4500, 4700, 4900; top: relative turbulent energy 
e-t/etotai] middle: relative turbulent pressure gradient 
Vp t |/|Vp tota/ |, {ptotai=P+Pt); and bottom: relative con- 
vective flux F c /F total . 



In order to perform the extensive survey of this pa- 
per we had to use a relatively coarse zoning of 180 mesh- 
points, viz. 40 uniform zones in from the surface, then 
increasing geometrically up to an inner temperature of 
2.5 million K. Just as for the radiative models, it is im- 
portant to provide a grid point right in the center of the 
sharp hydrogen PIR (unless a very fine mesh is used). We 
accomplish this by constraining the equilibrium model to 
have one anchor zone at (T = 11, 000 K), holding this an- 
chor fixed in the iteration to a full TC model. Failure to 
do so will result in erratic fluctuations of the growth-rates 
along sequences of models. We also find that for models 
with strong convection which penetrates deeply into the 
stellar envelope it becomes necessary to refine the zoning 
in the Fe PIR, i.e. around 200,000 K. 



We turn now to the properties of the equilibrium mod- 



els. 
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Fig. 2. Effect of a c = 
= 0.375, M=5M , L - 
et/etotai] middle: \Vp t \ 



5,4,3,2; Cepheid models with a\ 
-- 2090L©, and T e// =4700K; top: 
'|Vp tota /|; and bottom: F c /F tota i- 



3.1. The Convective Flux 

Along a Cepheid sequence (with fixed M, L and compo- 
sition) the size of the negative entropy gradients increases 
when T e ff decreases, causing convection to become more 
vigorous and to penetrate deeper into the stellar enve- 
lope. This is clearly illustrated in Fig. [l] for our stan- 
dard sequence of models with decreasing T e yy, for which 
we adopt the parameters M=5M Q , L=2090Lq, and with 
T e// =4300, 4500, 4700, 4900 K, respectively. Throughout 
this paper we use the compositional parameters X=0.70, 
Z=0.02. Similarly for the tests of this paper we adopt the 
standard TC parameters, viz. a c —h., at — I., a u =3., a\ 
=0.375, a s =0.75, a<j=4., and a p =2/3. These values are 
just used for reference purposes and are not meant to be 
the values to be adopted ultimately. A detailed calibra- 
tion of these parameters using data from observations is 
in progress and will be addressed in a subsequent publi- 
cation. 

In Fig. [I] (top) we display the fractional turbulent en- 
ergy e t /e to tai for the model sequence, where e to tai = e+e t . 
The behavior of the corresponding turbulent pressure ra- 
tio Pt/ptotai is similar, but slightly broadened because of 
the p dependence of pt- In the middle figure we show the 




Fig. 3. Effect of a A = 0.375, 0.328, 0.281 and 0.211; 
Cepheid models with a c =3, Af=5, L — 2090, and 
T e/ /=4700K; top: e t /e; Vpt|/|Vp tota z|; and bottom: 

F c j Ftotal- 



relative turbulent pressure gradient Vp f |/| Vp tota /|. Even 
though the turbulent pressure generally does not exceed 
a few percent, its gradient plays a more important role in 
the hydrostatic balance of the convection zones. The bot- 
tom figure displays the fraction of total heat flux that is 
carried by turbulent convection. In the hydrogen PIR it is 
close to 100%. 

In Fig. the temperature range 4300 - 4900 K corre- 
sponds roughly to the width of the IS for these models; the 
coolest model, near the red edge, is seen to have the most 
widespread turbulent pressure and the deepest and most 
efficient domain of convective flux. We point out that it is 
the least convective models that have the largest turbu- 
lent pressure gradients. This is an indirect result of those 
models being less turbulent and consequently less diffu- 
sive. 

The convection zones are of course associated with the 
PIRs. In Fig. |l| the H, He + , and He ++ convection zones 
are all visible. The maxima of the convective peaks stay 
approximately at the same temperatures, but the struc- 
tural changes of the models with decreasing T e ff causes 
them to move to higher pressures. For the cooler model, 
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the convection zone extends from the H PIR all the way 
to the Fe peak PIR. 

The convective flux F c (Eq. [?]), depends on the pa- 
rameter a c which we have held fixed in the sequence of 
Fig. § We now consider the changes that occur when a c 
is increased for a model with fixed T e ff = 4700 K and 
with the same M and L as in Fig. 0. Interestingly, Fig. | 
shows that a decrease in the strength of the convective 
flux (a c ) is remarkably similar to increasing the effective 
temperature of the model. 



3.2. Mixing Length 

The mixing length is the characteristic length-scale of the 
turbulent eddies. Therefore cha appears in a number of 
places, and in particular in the turbulent diffusivities, and 
ckA represents a measure of the nonlocality of turbulence. 
Its effect on the spatial distribution of turbulent pressure 
and convective flux is dramatic, and qualitatively differ- 
ent from that which results from changes to a c or T e ff. In 
Fig. U we see that the turbulent quantities energy (top), 
and pressure gradient (middle) become very localized to 
the unstable regions even when a\ is made only slightly 
smaller than the standard value of <xa — 0.375 that was 
adopted in Figs. Concomitantly, the region of efficient 
convective flux (Fig. ||, bottom) separates into increas- 
ingly narrow regions centered on the convectively unstable 
PIRs. 



3.3. Properties of Convection 

Efficiency - There are two conflicting definitions of con- 
vective efficiency. The Peclet number is one commonly 
used to measure the relative efficiency of convective over 
conductive transport, and is given by Pe = td/t c , or the 
ratio of the thermal diffusion timescale to the convective 
timcscale. Here, diffusion is radiative and its timescale is 
given by td ~ d 2 /x, where d — a^Hp is the mixing length 
and x — (16crT 3 ) / (3c p p 2 k) is the usual radiative diffusiv- 
ity The convective timescale is then given by t c ~ d/e\ /2 . 
The Peclet number for our models is large in the con- 
vective region (bottom Fig. ^), and convection is indeed 
efficient in the sense that it carries a large fraction of the 
total heat flux in the PIRs (Figs. 0- |). 

Truly efficient convection ought to remove the unsta- 
ble gradient that produces it. In Fig. [|, top and middle, 
we show the density and entropy profiles of two convective 
models. We see that in fact convection does not manage 
to remove the density inversion. It reduces the size of the 
negative entropy gradient by an order of magnitude com- 
pared to the radiative model, but the remaining entropy 
gradient is still strong, and it is shifted to the wings of the 
original convection zone. 



The Nusselt number Nu - This dimensionless number rep- 
resents the relative importance of convective and conduc- 
tive fluxes, and is defined as Nu = F c /F con d for a system 
that is ordinarily conductive, but here it is radiative. Nu 
is thus connected with the quantity plotted in Figs. _ |3[ 

The Nusselt number must be related to the Rayleigh 
number, which is the usual measure of the degree of con- 
vective instability, defined as Ra = g f)d?TY / (yx) where 
g = GM r /r 2 is the local gravity and v is the kinematic 
viscosity. Thus Nu ~ Ra a , and, under idealized condi- 
tions good theoretical arguments have been give n that , 
depending on the physical assumptions (e.g. Spiegel 1971 ), 
a = 1/3 or a = 1 /2, but for which no experimental support 
has yet been given. 

In Fig. ||we display the run of the local Nusselt number 
versus the local Rayleigh number (assuming Pr = 10~ 9 ; 
this value has essentially no effect on the slope) through- 
out two Cepheid envelopes; the lines represent the same 
models that appear in Fig. ||. The slope lies between 0.45 
and 0.53 throughout the convective region, a value that 
agrees well with the theoretical upper limit, considering 
the strong gradients and rapid variations of the physical 
quantities (of x among other things) in the convective re- 
gions. 



4. Linear Stability Analysis 

For the stability analysis we perform a linearization of all 
three Eq s. (|l| , ||, ||) , with that of all their related quantities 
Eqs. (p|-|lO|). As pointed out earlier, the equilibrium model 
depends on the a's (except a„), and consequently need to 
be reconstructed for any change of a parameters before 
the linear problem is solved. 

Because of the turbulent energy equation (Eq. ||) there 
is an additional spectral branch of eigenvalues compared 
to the vibrational and thermal spectra of the nonadiabatic 
radiative models. We compute the complete spectrum of 
eigenvalues by using canned eigenvalue solvers (e.g. RG 
or EIGRF, as suggested by Glasner & Buchler |1993| ). If 
only the lowest vibrational eigenvalues are desired it is of 
course much faster to first compute the adiabatic eigen- 
values (symmetric band matrix) and insert these as initial 
guesses into a Castor type iteration (Castor 1971 ). 

In this section we examine the role that the various 
ingredients of the equations play in the stability of the 
low- lying vibrational modes of Cepheid models. The dom- 
inant effects are the nonlocality (mixing length), and the 
strengths of the convective flux and of the eddy viscosity, 
as we will show. 

Cepheids are known to pulsate primarily in the funda- 
mental or first overtone modes, so we restrict our attention 
in the following analysis to just these two modes. The ef- 
fects of turbulent convection on the strange modes will be 
addressed in a companion paper. 
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Fig. 4. Structure of a typical Cepheid model (M=5Mq, 
£=2090L q ); soiid Jine: T e// =4900K; dotted line: 
7e//=5300K. The surface is on the right. All a's have 
standard values. Top: density (Log) ; dashed line: radia- 
tive model (T e ff= 5300 K); middle: dimensionless entropy 
gradient Y; bottom: convective efficiency (Peclet num- 
ber). 



4-1. Eigenf unctions 

Showing the effects of TC on all the lowest vibrational 
eigenvectors would take us too far. In Fig. ^| we limit our- 
selves to the modulus of the (complex) radial displace- 
ment eigenvector |5r/r| for two convective models (same 
as for Fig. ||). The purely radiative model (same M, L and 
T e ff= 5300 K) is shown as a dashed line. Convection re- 
moves some of the bump that occurs in the H PIR because 
of the density inversion (cf. Fig. ||) . 

The bottom figure shows the Lagrangean perturbation 
of the specific turbulent energy. As expected Set is local- 
ized to the convective regions. 



4-2. Work-integrands 

In Fig. [?] we show the work- integrands (oc Sptotai x Sv*) for 
two models (top, T e// =4900 and bottom, T e// =5300K) 
of Fig. ^j. Shown are the exact (or total) work-integrands 
(solid lines), together with the separate contributions of 



Log Nu 




Log Ra 



10 



12 



14 



Fig. 5. Nussclt versus Rayleigh numbers throughout the 
convective regions (where Nu > 1) of two typical Cepheid 
envelopes with M=5M , L=2090L Q , T e// =4900K {solid 
line), T e// =5300K {dotted line). 



the eddy pressure (long dashes) and the turbulent pressure 
(dotted line). The eddy pressure is always damping, as 
one can also show theoretically. The turbulent pressure 
can be both driving and damping, but it always plays a 
subdued role in the overall growth-rate, at least for the 
lowest frequency modes. 

The connected crosses show the decoupled case, i.e. the 
work-integrands with the variations of all the TC quanti- 
ties are omitted. Note that this is not equivalent to omit- 
ting both the Spt and 8p v , because the Sp includes a con- 
tribution from SF C and SC; cf. Eq. (|J). 

The short-dashed line shows the work-integrand when 
the variation SF C alone is omitted from the linearization. 



4-3. Mixing Length 

As we would expect from the behavior of the static models, 
the properties of the pulsational modes are also found to 
be sensitive to the value of the mixing length parameter 
a a ■ In Fig. || we show the growth-rates of the fundamental 
(top) and first overtone (bottom) modes in four sequences 
each with temperature ranging between 4000 - 7000 K, 
but having different a\. We define the relative growth- 
rates as those of the model energies, i.e. rjk = AnKk/uJk, 
obtained from the eigenvalues = iu>k + Kk- 

The mixing length is the fundamental length-scale in 
TC models such as this one. The value of oia therefore 
affects many of the model quantities that have to be kept 
in mind when interpreting the results of the variation of 
a a ■ Here we perform our study by varying only a a while 
holding all other a's fixed. 
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Fig. 6. Fundamental eigenfunctions for the two Cepheid 
models with M=5, L=2090, T eff = 4900 (solid) and 
5300 K (dotted) as in Fig. || Top: radial displacement 
vectors (\Sr\/r); dashed line: radiative model; bottom: 
turbulent energy eigenfunctions |<5et|. 



4900 K 


, i , , , , 

{ I'f 


-u,-,,,,.!!: 


\ K 

, , , , , , i , 


i j 

1 , , , , > 


/ 

/ 

1 

1 
i 
1 
1 

- 



5300 K 




exact -ssasi^ | 

P v part jf 

p t part \f) 

F c part 

- * ■ * * decoupled 

i , , , , i 


M. sSi^ \ \ I / 

\ I \* \ 
\' J i 
1 ' 
1 ' 
\ 1 
\l 

1 , , , , 1 



Log T 




7000 



T efl (K) 



Fig. 7. Work-integrand and its constituents for the mod- 
els of Fig. § Top: T e// =4900K and bottom: T e// =5300K. 



Fig. 8. Relative growth-rates for a Cepheid with fixed M, 
L and T e // in the range: 4100-7300 K. Top: fundamen- 
tal mode; bottom: first-overtone. Individual curves corre- 
spond to different values of q;a as indicated in the upper 
panel. 

The most unstable sequence corresponds to the small- 
est value of the mixing length (most local). 



4-4- Convective flux 

While the behavior of the cooler models and of the higher 
frequency modes is strongly tied to the mixing length the 
same is not true of the convective flux parameter ct c . In 
Fig. ^| we show again the fundamental and first overtone 
growth-rates for the same star, but now with «a = 0.375 
and a c = 2., 3., 4., 5.. 

We see here that while a\ has an enhanced effect on 
both the cooler models and the overtone mode, a c can 
affect the locations of both the red edge and the blue edge 
and the location of the maximum growth rate. 



4-5. Eddy Viscosity 

The eddy pressure represents the stress exerted on the 
fluid by the turbulent eddies - it can be characterized 
in terms of the eddy viscosity parameter. We can readily 
examine the effects of this term on the damping of pul- 
sational modes by varying the parameter ot v . In Fig. [l^ 
we show the fundamental and first overtone growth-rates 
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when a v = l.,2.,3.,4. (a v = 3. is what we take as the 
standard value). 

The effect on the overall shape of the IS is dramatic, 
and for this reason, each of the parameter dependences 
above (that of a c and that of cka) is presented for a range 
of eddy viscosities a u — 1., 2., 3., 4. . 



4-6. Limits: Adiabatic Turbulence and Instantaneous 
Turbulence 

Numerically we also have the ability to disallow perturba- 
tions to individual turbulent quantities; this is useful for 
determining which TC components, when perturbed, play 
the largest role in the growth-rates, for example. We now 
use this tool to examine some limits for the coupled con- 
vection and pulsation problem for our model with M=5, 
L=2060. Of particular interest are the adiabatic turbulent 
and the instantaneous turbulent limits. 

We can constrain the turbulent perturbations to be 
adiabatic in that the specific entropy of the turbulent ed- 
dies remains unchanged, i.e. Sst = 5e t — ptSv = 0. Note 
that this is not equivalent to ignoring the perturbation 
of Set, nor is it equivalent to assuming a very long TC 
timescale. With the help of Eq. || one sees that this limit 



allows for the finite compressibility of the turbulent en- 
ergy, e t ~ p Qp , or j t = 1 + otp). (We stress that our eigen- 
values are fully nonadiabatic though, as far as the usual 
gas and radiation entropies are concerned). These adia- 
batic turbulent growth-rates can be seen in the dotted 
curve of Fig. 11 compared to the full nonadiabatic tur- 
bulent growth-rates (solid line) ; the effect on these lowest 
modes is not negligible showing 'adiabatic turbulence' to 
be an unacceptable approximation. 

At the other extreme, we can allow the turbulent en- 
ergy perturbations to respond instantaneously, i.e. the 
turbulent timescale is taken to be much shorter than the 
pulsation time. In Eq. || this corresponds to setting the 
Lagrangean variation of the flux and source terms equal 
to zero (note that p t is part of the turbulent entropy, and 
p v disappears in the linearization). These instantaneous 
turbulent growth-rates, as seen in the bottommost curves 
of Fig. |n], lie very close to the exact ones. This indicates 
that overall, and at least for the fundamental mode and 
first overtone, the turbulent convection response-times re- 
main much shorter than the period of the pulsation. (If 
one wanted to simplify the linear stability analysis for the 
lowest modes in Cepheid models, one could thus drop the 
time-dependence in Eq. ^ and express Set in terms of the 
variations of the other quantities). 
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It is also interesting to ask which of the perturbed 
TC quantities contribute most significantly to the overall 
driving and damping of the pulsational modes. We have 
therefore added four more curves in Fig. 1 1 . 
1. The long-dashed curves show the growth-rates with the 
neglect of the eddy viscosity perturbation (5p„=0). With- 



out this critical stabilizing contribution the width of the IS 
would be excessively large. Comparison of this curve with 
those in Fig. |l0| clearly shows that as the eddy viscosity 
is reduced, the red edge in particular is rapidly shifted to 
unrealistically low temperatures. 

2. The short-dashed curves represent the growth-rates with 
the neglect of the Lagrangean convective flux perturbation 
(5F C =0). This term is seen to play an equally important 
stabilizing role. 

3. The connected crosses show the growth-rates corre- 
sponding to the neglect of the perturbations of all TC 
quantities. The effect of turbulent convection appears only 
in the structure of the static model. Not only are modes 
much less stable overall, in addition, the sequence does 
not have a red edge within a reasonable T e //- 

4. The large heavy dots represent a purely radiative se- 
quence, displaying a continuing instability for the coolest 
models (no red edge with a reasonable T e //). 




Fig. 12. The interpolated data surface, here for the fun- 
damental mode for M/M©=3, 4, 5, 6, 6.5, 7 (for M-L 
relation, cf text) from which each IS is determined. The 
t]q=0 contour superposed on the surface corresponds to 
the red (right) and blue (left) edges. 



Most of the period shift in the TC models as com- 
pared to their radiative equivalents is due to the structural 
rearrangement that convection causes in the equilibrium 
model. These structurally induced changes vary from less 
than 1% near the blue edge to as much as 15% near the 
red edge, for the very most turbulent models. 

The effects of the perturbations of the turbulent energy 
equation on the pulsation can also cause a shift in the 
linear periods (in a way that is similar to the nonadiabatic 
temperature perturbations). These TC shifts are generally 
less than one percent when compared to the nonadiabatic 
nonturbulcnt models. 



5. Instability Strips 

The shape and location of the IS depends on the a pa- 
rameters that ultimately we want to determine from the 
observational constraints. To this effect we now calculate 
Cepheid model sequences with masses in the range M=3 - 
7Mq. The composition is taken to be X = 0.70, Z = 0.02, 
and L is determined from an assumed mass-luminosity re- 
lation logL = 0.83187 + 3.56 log M. At each M (and L) 
a sequence is calculated with T e ff ranging from 4000 to 
7000 K in steps of 100-200 K. At all M and L(M) the 
T e ff dependent growth-rates functionally resemble those 
of Fig. || 

The global picture can be visualized in Fig. which 
shows a typical surface rjo — 7yo(L,T e ^) (obtained from a 
cubic spline interpolation of all the computed data) . The 
blue and red edges are shown as contour lines. 

We turn now to the morphology of the IS in theoretical 
HR diagrams as a function of some of the a parameters. 
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Fig. 13. Linear instability strips as a function of a v and a\ with fixed a c — 3.0, fundamental IS (light shaded area) 
and first-overtone IS (dark area). Bottom to top: a v = 1.5, 2.0, 2.5, 3.0; left to right: a A = 0.281, 0.328, 0.375, 0.390. 



Fig. |l3| shows the results for the mixing length pa- 
rameter aA (increasing rightward with values 0.288, 0.328, 
0.375, 0.390) and a v (increasing vertically with values 1.5, 
2.0, 2.5, 3.0), and with a c —3. (the a p , at, a s , a<j are kept 
at their standard values). The light-shaded regions form 
the fundamental IS (770 > 0) and the dark-shaded ones 
the first overtone IS (771 > 0). (The dots represent calcu- 
lated models, but some of the sharp features are a result 
of the interpolation and of the relatively crude numerical 
mesh). One notes that the sensitivity to a a is very large. 

Too weak a convective flux parameter (small a c ) gives 
too broad an IS, both for the fundamental (F) and for the 
first overtone (01). However, a small mixing length has 
a similar effect. On the other hand, too large a convec- 
tive parameter kills off 01. The figure shows that various 
combinations of parameter values provide approximatively 



equivalent IS's, e.g. running approximately from the top 
middle toward the bottom right. 



Fig. |14| shows the morphology of the IS in theoretical 
HR diagrams as a function of the convective flux parame- 
ter a c (increasing from left to right with values 2, 3, 4, 5) 
and a v (increasing from the bottom to the top with values 
1.0, 1.5, 2.0, 2.5), and with a A = 0.375 (the a v , a 8 , at, ad 
again are kept at their standard values). This figure again 
shows that various combinations of parameter values pro- 
vide approximatively equivalent IS's, e.g. running approx- 
imately from the left top toward the bottom right. 

The first overtone IS closes off at large L for some 
range of a values, which is in agreement with observation. 
Thus the overtone Cepheids with the longest observed pe- 
riods (Pi) give a strong constraint on the value of the a's 
(and on the M - L relation) . For reference we have drawn 



P.A. Yecko, Z. Kollath & J. R. Buchler: Turbulent Convective Cepheid Models: Linear Properties 



11 





Fig. 14. Linear instability strips as a function of a v and a c with fixed qa = 0.375, fundamental IS (light shaded 
area) and first-overtone IS (dark area). Bottom to top: a v — 1.33, 1.67, 2.0, 2.33; left to right: a c = 2, 3, 4, 5. 



iso-period lines Po = 5 and 10 d (solid lines) and Pi = 
2 and 5d (dotted) in Figs. [l3| and [l4|. Overall, changes 
in the values of the a c and a\ cause small shifts of the 
isoperiod lines. As already mentioned most of these shifts 
are due to the structural rearrangement caused by convec- 
tion. That is why a u , which does not appear in the static 
equilibrium equations, has almost no effect. The periods 
are only slightly affected by convection near the blue edge, 
but the pulsation-convection coupling can cause more ap- 
preciable changes near the red edge for the more turbulent 
models. 

We conclude that, just from a visual inspection of 
these IS diagrams, it is possible to delineate the acceptable 
ranges of the parameters a a, a v and a c . Roughly speak- 
ing the parameters with reasonable looking IS's lie on a 



2D surface in the (a c , a„, a\) parameter space. Table I 
quantifies this visual result over a more extensive range of 
the a parameters. 

This flexibility in the a values appears fortunate be- 
cause ultimately we want not just to fit observed funda- 
mental and first overtone IS's for the Galaxy, but also 
for the Magellanic Clouds which have lower Z values. In 
addition, the wealth of observational Cepheid Fourier de- 
composition parameters of the (nonlinear) light and radial 
velocity curves in these galaxies will need to be fitted. 

We have thus narrowed down the parameter space on 
which more detailed calibrations will be made that take 
into account the wealth of observations. However, we have 
not made any IS diagrams showing the effects of the other 
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Table 1. Examples of three parameter combinations 
which give nearly equivalent results for the F and first 
overtone instability strips. 





Cth 


Ct v 


0.5 


0.375 


1.50 


0.5 


0.500 


1.00 


3.0 


0.150 


9.00 


3.0 


0.328 


3.00 


3.0 


0.420 


1.00 


10.0 


0.100 


25.00 


10.0 


0.200 


6.00 



four a's because they seem to play a lesser role, at least at 
the level of accuracy considered in this parameter survey. 

6. Conclusions 

Purely radiative Cepheid models cannot satisfactorily ac- 
count for all observational data, and this has forced us to 
include turbulence and convection in our Cepheid models. 
We have adopted a simple ID model diffusion equation for 
turbulent convection (TC) that can easily be incorporated 
into a radial stellar pulsation code. Although there is no 
guarantee that such a single ID equation can satisfactorily 
approximate the effects of turbulence and convection it is 
worth exploring to what extent we can obtain agreement 
with observation. The model TC equations contain several 
dimcnsionlcss, order unity parameters that are directly re- 
lated to the physical quantities of the model that need to 
be calibrated with the help of observational astronomical 
constraints. 

Instead of performing time-consuming fully nonlinear 
hydrodynamic computations of the pulsations, we have 
developed a code to perform an efficient linear nonadia- 
batic stability analysis of the models. This has allowed us 
to make an exploration of the sensitivity to the free pa- 
rameters of the TC equation, both of the structure of the 
equilibrium Cepheid models and of their linear properties. 

In agreement with other studies we find that the static 
Cepheid models exhibit convection primarily in their H 
and He partial ionization regions (PIRs), although in some 
parameter range a convective zone also appears in the Fe 
PIR. For strongly convective models the separate convec- 
tive zones merge into a large one that penetrates up to 
320,000 K. 

The coupling of pulsation and convection also modi- 
fies the linear properties of the pulsational modes, as ex- 
pected. It also creates a new branch of turbulent diffu- 
sion modes linear that our stability analysis establishes to 
be extremely damped, fortunately for our ID TC recipe. 
These additional modes thus cannot create havoc numer- 
ically or physically. This large damping is consistent with 
the convective timescales being very short compared to 



the periods of the vibrational of the excited modes which 
themselves are short compared to the growth-times. 

Our survey shows that the stability of the fundamental 
and first overtone modes depends on all the TC parame- 
ters, but that it is dominated by (1) the convective flux, 
(2) the eddy viscosity, and (3) the mixing length. Gener- 
ally, the effects of eddy viscosity and mixing length are 
greater in the higher overtones. 

In agreement with other much less detailed TC studies 
we find that reasonable a parameter values exist for which 
IS's of the galactic Cepheid F and 01 have an acceptable 
range of temperatures. In this paper we have made no 
effort to narrow down the acceptable range of TC param- 
eters by using the large amount of observational informa- 
tion. An application of the TC code to a broader sample of 
Cepheid models and its comparison to the large data base 
of Galactic, and Magellanic Cloud Cepheids is in progress, 
as is an extension to RR Lyrae stars. 
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